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Abstract. In this article a continuous time deterministic model with vaccination and treatment strategies is 
formulated to assess the effect of reinfection on the transmission dynamics of Tuberculosis (TB). The 
involvement of reinfection in our model causes relapse and leads to the possibility of backward bifurcation at 


critical value of effective reproduction number R, = land hence the existence of multiple equilibria when 
effective reproduction number R, <1. This indicates that even by reducing effective reproduction number 
R, below one is no longer a sufficient condition to eradicate the disease from community. An additional 


reduction of effective reproduction number R, below the saddle-node bifurcation value is required to 


eradicate disease from community provided that the disease free equilibrium is globally asymptotically stable. 
Numerical simulation results are presented to validate analytical results. We suggest that reinfection is an 
important feature of TB and has to be considered when modeling the complex dynamics of TB. 


1. INTRODUCTION 

Tuberculosis (TB) is a chronic bacterial infectious disease caused by pathogen 
Mycobacterium tuberculosis with more than one-third of the world human population as 
its reservoir [1, 9, 16]. A global annual estimate of 8.6 million people develop Tuberculosis, 
of which 1.3 million die from disease. It is reported in [24] that, the burden of disease 
caused by TB is high in developing world where poor nutrition, congested accommodation 
and emergency of HIV are manifested. The global estimates of incidence, prevalence and 
mortality rates per 100,000 population in 2012 were respectively 255, 303 and 26 and 
Tanzania incidence, prevalence and mortality rates per 100,000 population were 165, 176 
and 13 respectively as per [24]. It therefore raises a quest to find desirable means to curtail 
TB morbidity and mortality rates. 

Tuberculosis disease is mainly of two types: pulmonary and extra-pulmonary TB. 
Pulmonary TB is a common form of TB that affects lung while extra-pulmonary TB affects 
other parts of body and organs including central nervous system and bone [23]. This 
particular study focuses on pulmonary TB. Tuberculosis is an epidemic disease spreading 
in the air when the infectious person with pulmonary TB expel bacteria by coughing, 
singing, sneezing, speaking and so on [6]. An individual with active TB has usual symptoms 
which are general weakness or tiredness, fever, weight loss, loss of appetite and night 
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sweats. Further symptoms are coughing, coughing up of sputum and/ or blood, shortness 
of breath and chest pains if the infection in the lung get worse [7]. TB draws back 
economics of the world and Tanzania in particular as it affects men than women and 
especially the productive working group [23]. In absence of HIV a small proportion of 
about 10% of infected individuals with Mycobacterium tuberculosis develop TB and 
becomes infectious within two years upon infected [20]. Most become latent for the rest of 
their lives as long as their immune system is not compromised [6]. The recovered 
individuals from TB do not acquire the permanent immunity. Some of them they become 
latent again. Even with treatment interventions, the rates of reinfection TB are higher than 
those of new TB [20]. Mathematical modeling of epidemiology of Tuberculosis has recently 
become the powerful tool to study the complex dynamics of the disease and explore the 
role of various TB features such as reinfection and reactivation. [9] formulated 
mathematical model of TB with exogenous re-infection. The results of their work suggest 
that exogenous reinfection has drastic effect on qualitative dynamics of TB and it allows the 
possibility of subcritical bifurcation at critical value of basic reproduction number œR, =1. 


Buonomo and Lacitignola [3] applied bifurcation method introduced in Castillo-Chavez and 
Song [6] based on the use of center manifold theory [4] to derive conditions for existence of 
either forward or backward bifurcation for vaccination model introduced in [10]. They 
clearly explain the role of vaccination, treatment and transmission parameters for the 
occurrence of forward or backward bifurcation. Okuonghae and Aihie [19] examined the 
effect of Direct Observation Therapy Strategy (DOTS) on dynamics of TB against the 
fraction of active cases detected. They formulated mathematical model that involves the 
fraction of detected cases undergoing treatment under DOTS and other fraction not 
detected. The qualitative analysis of this model shows that in presence of exogenous re- 
infection, reproduction number must be outside the bifurcation range for disease free 
equilibrium to be asymptotically stable. The parameter for case detection has shown to be 
important in reducing backward bifurcation range as well as reducing the reproduction 
number. They further argued that if the critical level for case detection parameter is not 
reached then TB persist in population and become endemic. [12] propose mathematical 
model of TB that includes exogenous reinfection in order to understand the recent increase 
of TB incidences in Korea. In their study, parameter for case finding effort was found to be 
significant impacting component on curbing active TB cases. They recommended that for 
dramatic reduction of TB incidences, the treatment of active TB cases should be 
accompanied by case finding (taking medication before the actual active TB has clinically 
diagnosed) effort than to take each measure alone. This paper concentrates on 
investigating the role of re-infection on the one strain TB model with vaccination and 
treatment and its impact on TB transmission dynamics. Bifurcation and stability analysis of 
equilibrium points are properly investigated. 
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2. MODEL FORMULATION 
Our population model is subdivided into six compartments and is developed from the basic 
SEIT (Susceptible-Exposed-Infectious-Treated) compartmental model. A compartment of 
Vaccinated population (V ) is added to form SVEIT model. In addition compartment of 
infectious population (Z ) is subdivided into two compartments which are severely infected 
population (/,) and mildly infected population (Z7, ). Severely infected population (/, ) 
progresses faster to treatment group compared to mild infected population (Z, ). In this 
model susceptible population will be recruited at a rate 4. Some susceptible individuals 
will come into contact with infectious individuals and being infected at a rate of B. A 


proportion, p of babies will be vaccinated at birth while the remaining proportion (1-p) 


will be left out of vaccination to join the susceptible population. Once vaccinated babies 
loose immunity they become susceptible at per-capita rate 0 , whereby 1/0 is the period 
after which a vaccinated baby looses immunity. The Latently infected individuals progress 
to active TB through endogenous reactivation. The proportion (l—7) of Latently infected 


individuals progresses fast to severely infected class, Z, while the remaining proportion, 7 
progresses slowly to mildly infected class, J, at the same per-capita rate e. Under usual 


circumstances mildly infected individuals take a long time to progress to treatment group, 
T than severely infected individuals. That is a proportion, b of mildly infected individuals 


progresses to treatment group, T while the remaining proportion, (1—) progresses to 


severely infected class, /, at the same per-capita rate w. The severely infected individuals 


progress to treatment group at a rate of v. The treatment group, T is assumed to undergo 
exogenous re-infection and relapse back to Latent group with infection level, y. The 


infectious individuals /, and Z, are assumed to die at disease induced mortality rates of ô, 
and ô, respectively while the rest die naturally at a rate of uw. All variables and parameters 


are assumed to be non-negative. 
In addition the following assumptions are taken into consideration during the formulation 
of the model: 
i. All individuals are born susceptible. 

ii. ©The members of population mix homogeneously. 

iii. Age, sex, social status, do not affect the probability of being infected. 

iv. Natural recovery is negligible and hence ignored. 

v. | Vaccinated population looses immunity and become Susceptible. 

vi. No more Vaccination can be administered to an individual infected with TB 


or to someone who previously was vaccinated. 
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vii. Once recovered from Treatment an individual reverts to be Latent and may 
experience another episode of disease. 
viii. Once an individual is infected he/she will not recover if no treatment is given. 
The above description of model formulation together with the assumptions leads to 
compartmental diagram in Figure 1.The full description of variables and parameters used 
to formulate the model are in Table 1 and Table 2 respectively: 


Table 1: Description of variables of the model. 
Variable Description 


S(t) The Susceptible who are at risk of being infected at time t. 
Lit) The latently infected individuals at time t. 
V(t) Vaccinated individuals at time t. 
I(t) Individuals who are severely infected with TB at time t. 
I, (t) Individuals who are mildly infected with TB at time t. 
T(t) Individuals Treated against TB at time t. 
Cee + 64 dy 


— 
l 1 
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Figure 1: Schematic flow diagram showing dynamics of tuberculosis, where I =Z +Z. 
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Table 2: Description of Parameters of the model. 
Parameter Description 


A Per capita birth rate. 

p Per capita infection rate. 

p Proportional of babies who are being vaccinated at birth. 

0 The rate at which a vaccinated individual looses immunity. 

E The rate of progression from Latent class to both severely and mildly Infected 
classes. 

1] Proportional of Latently infected population progressing to mild infected class. 

H Per capita natural death rate. 

Ô, Per capita additional death rate of severely infected class. 

Ô, Per capita additional death rate of mildly infected class. 

Q Proportional of mildly infected class who are treated. 

w The rate at which a mildly infected individual is transferred to both severely 
infected and treatment classes. 

D The rate at which a severely infected candidate is transferred to treatment 
class. 

y The factor that reduces the level of reinfection. 

2.1 Equations of the Model. 


Basing on assumptions made and relationship that exists between variables and 


parameters shown in Figure 1 the system of six ordinary differential equations that 
describes the dynamics of tuberculosis in presence of vaccination and treatment is given by: 





dS (1,+1,) 

Bap) N pS ag ey 

(l-p)N-£ — 

dV 

A Vea 

PN -(H+8) 

dL (+L) (+Z) 

oe (u+e) 

T =(1-4)eL+ (1-¢)ol, = (4+8, +091, (1) 


To =yeL-(u+0+8,)I 


dT (+L) 
= aol dhol —| uri =T 
7 i+ ool, c yp x 


2 


N=S+V+L4+if+1,+T. 


By adding the state equations in (1) we end up with rate of change of population, 
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= (4-n)N-8, 5, (2) 
2.2 Normalization of the model. 

The model (1) can easily be analyzed after being normalized such that the total population 
is one. The normalization is done by scaling the population of each compartment by total 
population. We transform the actual proportions by setting: 


S V, L, L.L 


E E E a (3) 
N N N N N N 
where by s+v+l+i +i, +h=1. 
Substituting (3) into (2) we end up with: 
T =(4-4-8i-8i,)N (4) 
í 


Upon differentiating the proportions in (3) with respect to time t and make simplification, 


leads to the following dimensionless system: 


Z (1-p)4+øv-(4+8(i +i,)-4i,-6,i,)s, 

Z = ph-(1+0-di -ôi )v, 

LN +i, )+yph(i, tL) (Arei -ĝi )l, 

dt (5) 
S =(1=n)el +(1-4)oi, (4+8 40-04, -8i )i- 

= el (140046, -òi =o l)i 


< = vi, + oi, -(A+ yB (i, +i, )- ĝi - ôi, )h. 
subject to condition s+v+/+i, +i,+h=1. It can be shown that all the feasible solutions of 
system (5) enter the region of biological interest defined by 
Q= {(s,v,L,i,i,h) ER’ :st+v4+ltit+i,+h= 1} 


that is positive-invariant. It is enough to consider the dynamics of the flow generated by 
system (5) in Q. In this region, the model (5) is considered to be both biologically and 
mathematically well posed [11]. 
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3. ANALYSIS OF A MODEL 
We analyze model (5) in order to get some insights on dynamics of TB disease and 


transmission. 


3.1 Existence and Local orm of DFE 
Let E, =(s°,v’,/',i,,,i,",4°) be a DFE point of model (5). We set zero to the right hand side of 


Bos L, > 
each equation in (5) and assume that in absence of disease attack, / =i, =i, =h=0 to solve 


the steady state solution. The disease free equilibrium point is therefore given by : 


„a _(0+40-p) pà re 
E; =(s°, vf ,i ,t, ,h )= — „0,0,0,0 |. Before we prove for local stability of 
172 A+0 (A+04° 


DFE we define and determine the effective reproduction number, R, of model (5). 
Definition 1. The effective reproduction number, R, is defined as the measure of average 
number of infections caused by a single infectious individual introduced in a community in 
Which intervention strategies (in our case is treatment and vaccination) are administered 
[18]. 

The effective reproduction number R is computed by using next generation operator 


method [22] and found to be: 
— B(0+A(1-p)) (1-7)(A+@+0, )e+(1-) ane NE 


= (040) «| Caesone esos) a) 


E (A+0)(A+e)(A+@+6,)(A+6, +0) 


(6) 


Theorem 3.1. 7he disease free equilibrium of model (5), given the effective reproduction 


number, R, is locally asymptotically stable if R <1 and unstable if R >1. 


We prove Theorem 3.1 for local stability of DFE by asserting that the trace and determinant 


of Jacobian matrix at DFE denoted by J (E o) are strictly negative and positive respectively. 


Jacobian matrix evaluated at disease free equilibrium point is given by: 
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-A Æ 0 (praia) (pra | 2HE) 


4+0 A+0 
0 -(A+@) 0 aun E 0 
A+0 A+0 

J(E,)= = = 

A+0 A+0 
0 0 (l-n)e —(A+6, +0) (1-do)@ 0 
0 0 NE 0 —~(A+@+0,) 0 
0 0 0 V da —À 


Trace and determinant of matrix J(E,) denoted by Tr(J(£,)) and det(J(E,)) are 


respectively given by: 
Tr(J(E,))=-(64+0+e+0+@+6, +ô,)<0 


and 
0+4(1-p) 0+A(1-p) 
ree) s| ize) s| ie) 
det( J (E,))=-# (4+0)(1-n)e  -(2+8, +0) (1-o)a 


NE 0 —(A+@+6,) 


B(0+A(1-p)) 
=—1' (4+8) (A+0) 
—(A+e)(A+@+6,)(A+6, +0) 


E pepe a E 


(4+8) (Ate)(A+a+06,)(A+6, +0) 


{(l-n)(A+@+6, e+((I-6)o+4+6, +0) ne! 


=—A(R, -1). 
where A=7'(A+0)(A+e)(A+6, +v)(A+@+6,)>0. 
We find that Tr(J (E,)) is strictly negative and det (J (E,)) is strictly positive if and only if 














R, <1. We therefore conclude that DFE is locally asymptotically stable. 


3.2 Global Analysis of DFE of a model with interventions 
We analyze the global stability of disease free equilibrium point of model (5) by using an 
approach presented in [5]. The model (5) can be written in the following format: 
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"= A(X,-X, ,)+AX,. 
(7) 





From (7), X, and X, are vectors of non-transmitting and transmitting compartments 
respectively. X, „ is a vector at disease free equilibrium point E, of the same length as X,. 


From model (5) we define: 


X, =(s,v,h) ,X,=(Li,i,) .X,, [are me | ‘and 


A+O 4+0 
0+ A(1 — p) 
A+0 
X,-X,,> ye . For global stability of DFE we need to show that matrix A 
| + 
h 


has real negative eigenvalues and A, is a Metzler matrix (i.e. the off-diagonal elements of 
A, are non-negative, symbolically denoted by A, (x,)2 0,Vi + j). Using system (5), then 


the first and second equations in (7) can be written respectively in expanded form as: 


0+A(1-p) 
Foie cai a 
(1- p)A+Ov-(A+ B(i, +i, )-6i, -8,i )s A+0 / 
pa—(A+0-6i, —6,i, )v =A v- +Ali 
vi, + boi, —(A+B(i, +i,)- i, -ôi )h - i, 
and 
ps(i +i, )+yBh(i, +i,)-(4+e-8i -ôi )l / 
(1-n)el+ (1-6) qi, —(A+6, +0-6i, -ôi )i, |= 4,] i, |. For compatibility, matrices A, A 
nel —(A+@+6, -ôi -ôi )i, L 


and A, should be of order 3x3. By using non-transmitting elements from Jacobian matrix 
of system (5) and representation in (7) we find that: 

—A 0 0 0 (-£+6,)s (—B+6,)s 
A=| 0 -(A+0) 0/|,4=|0 OV OV and, 


0 oO -J 0 v+h(d,-yB) ba+h(d,-78) 
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—(A+e) B(s+yh)+o6]l B(s+yh)+6,l 
A, =| (I-n)e —(A+0+6(1-i)) (1-))@+6,i, |. We find that A is upper triangular 
NE Ot -(A+@+6, (1-i,)) 


12 
matrix whose eigeinvalues are located on its main diagonal. Therefore eigenvalues of A (i.e. 
—-4,-(4+0) and —A) are real and negative. In addition A, is a Metzler matrix since its off- 
diagonal elements are non-negative. That is 0</,i,,i,,n<1 and both (1-i,) and (1-i,) are 


strictly positive. Therefore DFE for system (5) is globally asymptotically stable in region Q. 
We have established important theorem: 

Theorem 3.2. The disease-free equilibrium point is globally asymptotically stable in Q if 
R, <1 and unstable if R >1. 


3.3 Existence of Endemic Equilibrium Point (EEP) of model with interventions 

Let E,(s',v’,l',i,,i,",4°) be an endemic equilibrium point of model (5). The conditions of 
existence of endemic equilibrium point E, are obtained by setting the right hand side of 
each equation in (5) equal to zero and solve model (5) in terms of force of infection 
f = (i, +i,) at steady state. Let k =A—6ji, —6,i, >0, for any pairwise choice of i," and i, 
values at endemic equilibrium. An endemic equilibrium point in terms of force of infection 
is given by: 


me a ee- r a eres ee 
Fae Tr tern Pees +k)(6,+0+k)—eyf" (va, +a,)} 
(F +k)(O+k){(e+k)(o+d, +k) (of? +k)(5, +0+k)—ayf" (va, +a,)} 
we define, a, =(1-7)(@+6,+k)+(1-6)on; a, = ban(d,+o0+k) 


l, 


(8) 
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If we substitute representations of i, and i, from (8) into the force of infection, 
f =B(i, +i; ) or f -B(i, +i, )=0 we find that: 
he(O+k(1-p))(of +k) f'| {(1-)(@+d, +) + (1-o) on} + (6, +0 +k) | 


= =0 (9) 
(F° +k)(O+k)\(e+k)(o+, +k)(yf"+k)(ô, +v +k)-ezf" (va, +a,)} 
Manipulating and simplifying (9) we end up with the following cubic polynomial: 
f (Af? +B f’ +C,)=0 (10) 


where by 
A, =| M-(vP+Q)(0+k)], 
B =k(M +y| M -(vP+0)(0 +k) |)-Zay(P+(ô, +o+k)ne), (11) 
C = k (Mk — fa, | P+(6, +o+k)ne |). 

Furthermore in terms of parameters of model (5) we define: 

M =(0+k)(e+k)(@+ô, +k)(6, +ot+k); P=(1-7)(@+6, +k)e+(1-o)ane; 

Q = bane(6,+0+k) and a, = (0 + k(1—p)) . We write C, in the following format: 


C,=(0+k)(e+k)(@+6,+k)(6,+0+k)k (1-R, ) (12) 


From (12), (0+k)(e+k)(@+06,+k)(6,+0+k)k* >0 and R, is effective reproduction 
number as indicated in (6). 

From (10), f° = (i, +i,")=0 corresponds to Disease Free Equilibrium (DFE) that we have 
already discussed while A f” +B f* +C, =0, that can be also be written in the form: 


2A 


1 


(13) 


Satisfies Endemic Equilibrium. The value of A is strictly positive. Depending on the signs 
of B and C, we have three cases to consider in order to have positive root of force of 


infection as follows: 

Case 1: In absence of re-infection we find that the parameter for level of reinfection, y=0. 
This implies from (11) that A =0. The polynomial A f” +B f* +C, =0 becomes linear, i.e. 
Bf'+C,=0or f = i, If B >0 then system (5) has stable endemic equilibrium when 


1 


C <0. This equilibrium happens when R, >las interpreted from (12). In this case 


backward bifurcation is not possible due to absence of multiple equilibria. 
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Case 2: Exactly one endemic equilibrium point. From (13), suppose B, <0 and C =O or 
B? —4AC, =0. This means the polynomial has just one positive root and hence the system 


(5) has unique endemic equilibrium. 
Case 3: Two endemic equilibria 
If B <0, C >0 and B°-4AC,>0, then the polynomial Af? +B f*+C =0 has two 
positive roots. This means that the system (5) has two endemic equilibria and hence the 
possibility of backward bifurcation. These three cases are summarized under the following 
theorem: 
Theorem 3.3: The number of positive endemic equilibria of Tuberculosis model (5) is 
hereunder summarized as follows: 

i. fC, <0, R, >1, the system has a unique endemic equilibrium. 


ii Uf B<O and C =0 or B’-4AC, =0, the system has exactly one endemic 
equilibrium. 

iii. Zf B <0, C, >0 and B? -—4AC, >0, the system has exactly two endemic equilibria. 

iv. Otherwise there are no endemic equilibria, Le. when AC, >0 and B >Q. 
From (iii), the critical point of effective reproduction number R° at which a backward 
bifurcation occurs is computed by setting the discriminant in (13) equals to zero. Thus, 
B? —4A C, =0 implies that 
B? -4A (O+k)(e+k)(@+6, +k) (6, +0+k)k’ (1-R*)=0 and, 

2 


C 


. Thus backward bifurcation occurs in the 


range R< R <1. Furthermore, we note from (13) that disease will be endemic if force 
of infection is strictly positive (i.e. f° >0) and both B and AC, are strictly negative. Thus, 
A >0 and AC, =A (0+k)(e+k)(w+6,+k)(6,+v0+k)k’ (1-R,)<0 if and only if R >1. 


* 


1 ? 


Therefore endemic equilibrium point E, (aa i, ,h’ ) is stable if and only if R >1. 


3.4 Stability of Endemic Equilibrium Point (EEP) of model with intervention 

The stability of an endemic equilibrium FE, of model (5) is analyzed by using Centre 
Manifold theory [4] as described in Theorem 4.1 of Castillo-Chavez and Song [6]. We 
change the variables of model (5) by setting s=x,v=x,,l=x,,i, =x,,i, =x,,h=x, such 
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6 
that 2% =1. We define vector X =(x,.%,,%).%,.%,%,) and F=(f. f, fo fo fo FY in 


dX 
such a way that the model (5) is re-written in the form T = F as follows: 


KE =(1-p)4+8x, -(4+ B(x, +x,)- ôx, -6,x,) x 
R =pA—(A+0-6,x, —ô,x, ) x, 

x, = f, = Px (x, +x, )+ ypx (x, +x, )- (4 +e- ôx, -—6,x, ) X, 
=, =(1-n)ex, +(1-6) ox, -(4 +ô, +0-6,x, —6,x,) X,, 


x, = f; = Nex, -(4+@ +ô,- ôx, —6,X,) Xs, 


P (14) 
X, = f, = 0x, + pox, -(A +B (x, + x, )-4,x, —6,%5) Xe- 
The Jacobian matrix J(E,) of system (14) at disease free equilibrium Æ, presented in 


Section 3.1 is given by 


—A 0 0 (-B+6)r (-Bt+6,)r 0 
0 - © + 0) 0 ON, Ôr, 0 
0 0 —(A+eé) pr pr 0 
J(E,)= 15 
S e Ime -(44640) ao 0 a 
0 0 NE 0 —~(A+@+0,) 0 
0 0 0 V A, —À 
0+4(1-p) pa l 
From (15) we define r =———~—— and r, =-~. In particular case when basic 
4+0 A+0 


reproduction number R =1, we choose our bifurcation parameter be p and consider our 

bifurcation to take place at f= f". Solving p from (6) when R, =1 we find that: 
(A+ 0)(A+e)(A+@+06, )(A+6, +0) 

(0+A(1-p))| (I-n)(A+@+6, )e+((I-d)o+1+6, +o) ne | 


The Jacobian of transformed system (14) at p =f" has simple zero eigenvalue that allows 


pap = (16) 


us to study the dynamics of the system (5) at p =f" using Centre Manifold theory [4]. The 
Jacobian of (14) denoted by J(£,) at 6 =f" has right eigenvector that corresponds with 


. z T 
zero eigenvalue given by œ = (w, wW, W, , W3, Wz, w, ) , where by: 
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É - B'rne)(àr O +(-B" +6,)(A+9)1,) + Brine (6,70 + (-B" +ô,)(4+0)r) 
Sic ee ee ee eee ee 


A(A+0) Brine _ 
ôr, ((4+e)(4+@ +ô, )- B'rne)+ inf rine 
W, m —— — —<—<—_-— __ wa... A- a >...) W; > 0, 
(A+0) Bre 
(A+@-+,) 
W, a * > 0, (17) 
m =| eiel Pit p T 
P rye 


w, = w, > Q, free. 


ae feh ee 
w =| Xo ET |, >. 
AP KNE 
, 6+A(1-p) 
From (17), k =(A+e)(A+@+06,)— frye and r Se By using (16) we show that 
+ 


k, is strictly positive justifying that the components w,,w,,w, >0 as follows: 


k =(4+e)(4+@+òô,)- prne 


-aro arosa) i E] 


=(At+e)(A+@+ o 0.0 | CA, ie 1 S 
=(At+e)(A sy SS] 
=(A+e)\(A+@t+ Cn a U 3 

=(/ \U at ae r 0. 


Moreover, the Jacobian matrix J (E, )at p = f" has left eigenvector ¥ = IAAL A AAN 
associated with zero eigenvalue satisfying the relation ¥ -œ =1, where by: 


COREE EN ee y _(li-o)o+4+6, +0) 
Br —— A+o+o, 


1 


Y =Y, =Y =0,¥, = Y, >0. (18) 


We compute the value of a and b that will govern totally the local dynamics of system (14) 
and determine whether it exhibits forward or backward bifurcation by employing Theorem 
4.1 of Castillo-Chavez and Song [6] and as restated in Theorem 3.4. 
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Theorem 3.4. Consider the general system of ordinary differential equations (14) with a 

parameter p such that a = Laps TR"x ROK and fe C’ (R' xR), where Q is an 
í 

equilibrium point of the system (i.e. f (0,8)=0 forall p ) and 


of. 
of 


j 


(0, 0) is Jacobian (linearization) matrix of the system around 


1. =p soo| 


the equilibrium 0 with p evaluated at O, 


2. Zero is a simple eigenvalue of A and other eigenvalues of A have negative real 





parts; 
3. Matrix A has a right eigenvector œ and a left eigenvector ¥ corresponding to Zero 
eigenvalue. 
th O° T k 
Let f, be the k" component of f and a = eww mm 0) and 


k,i, j=1 Ox.O 





b= yvw ,0) then the local dynamics of the system around the equilibrium 
k ,i=1 

point O Is mn determined by the signs ofa and b. In particular, ifa>Q and 

b>0O then a backward bifurcation occur at p =0. Signs of a and b play the vital 

role in describing the local dynamics of model (14) around equilibrium point Q as 

follows: 


a) a>0,b>0, when <0 with |$|«1, 0 is locally asymptotically stable and there 
exists a positive unstable equilibrium, when 0<\|f|«1,0 is unstable and there 


exists a negative and locally asymptotically stable equilibrium. 
b) a<0,b<0, when B<0 with |B|«1, 0 is unstable when 0<|f|«1,0 is 


asymptotically stable and there exists a positive unstable equilibrium. 
c) a>0,b<0, when B <0 with |f|«1,0 is unstable and there exists a locally 


asymptotically stable negative equilibrium, when 0<|f|«1,0 is stable and a 


positive unstable equilibrium appears. 

d) a<0,b>0 when f changes from negative to positive, O changes its stability from 
stable to unstable. Correspondingly, negative unstable equilibrium becomes positive 
and locally asymptotically stable. 

Computation of a and b 
We compute the value of a and b that will govern totally the local dynamics of system (14) 
and determine the conditions for existence of backward bifurcation following the signs of 
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a and b by employing Theorem 4.1 of Castillo-Chavez and Song [6] and as implied in 
Theorem 3.4 of this article. 


Since the components of left eigenvector ¥ =¥% =¥ =0 (for k=1,2 and 6) we compute 


the values of a and b for only k =3,4,5. The only non-zero second order partial derivatives 
of (14) at DFE when f= p° are: 

















Of, Of, =p Of, _ Of, eg OSs 895, AA 
Ox Ox, 0x, 0X, Ox,0X,  0x,0X; Ox, Ox, 
n Of. 
By using 4= ) P ww, (0,0) we compute a as follows: 
ae OX,OX, 


a=V¥, (w,w,8 +ww.B +Ww,w,0, +w,w.0, +w,w yp +ww yp ) 


+, (2w,’6, + w,w,6, ) + (w,W,6, +2w Ò, ) (19) 
— Pp w, (w, TU ) a wy Pp (w, TW; ) + W,0, (Pw, ij 2P,w, E Pw, ) 
+ w.6, (Pw, + ¥,w, + 2%, ). 
| , os, 
On the other hand, the value of b is computed by using the formula, b= > Pw, (0,0). 





k i=l i ax. op” 
The associated non-zero second order partial derivatives of (14) at DFE when f= f and 
k =3,4,5are: 


OF 0 Ise OP) 








f EA The value of b is therefore given by: 
oxp  Ox.0p A+0 
0+A(1- 0+A(1- 6+A(1- 
a ,| were) PA) yf PHA) sw)o0 (20) 
A+0 A+0 A+0 


From (19), let ¢, = ¥,8"w, (w, +w, ) and 
C, =w y PP (w, +w, )+ wò, (Ew, +28, w, + Pw, )+ wô, (Pw, +¥,w, +2¥w,). It follows that 
the sign of a depends on the value of w. If w >0 or w <0 and ¢, >¢, then a>0. We 


formulate the following theorem. 
Theorem 3.5: w >0 or w <0 and Ẹ, >¢,, a>0 then model (5) exhibits backward 


bifurcation at R =1. If f <0 then there exists unstable positive endemic equilibrium point 
and correspondingly if p >Q then there exists a stable negative endemic equilibrium point. 


Therefore endemic equilibrium point is locally asymptotically stable if R >1 but close to 1. 
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BAChKVYARD BIFURCATION. 


Stable EE 


Stable DFE 


Unstable DFE 


Force of infection 


Unstable EE 





o 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 
Effective reproduction number, R_ 


Figure 2: Bifurcation diagram showing backward bifurcation with estimated parameters 
B=14 5 y=1.8; 0 =0.8 ; €=0.396 ; n=0.1;A=0.9; 5 =0.3; o=0.6; p =0.1; v=0.9;36, = 0.2 
w=0.6; and o=0.1 for numerical simulation. 

Figure 2 shows the backward bifurcation of system (5) that occurs at threshold parameter 
R, =1, due to presence of multiple equilibria and re-infection. DFE stands for Disease Free 
equilibrium and EE stands for Endemic Equilibrium. In the neighborhood of 1 when R <1 
then stable DFE coexists with two endemic equilibria: the small unstable EE (with smaller 
number of TB infectives) and larger stable endemic equilibrium with large number of 
infectives. This implies that even with classically reducing the threshold parameter to less 
than unity does not clear TB from community. That is why we say backward bifurcation is 
an undesirable feature of TB. When R, >1 then we have two equilibria: unstable DFE and 
large stable EE. According to Buonomo and Lacitignola [3] if R, is nearly below one then 


disease control depends on initial sub-populations of the model under consideration. That 
is reducing R, below the critical value R‘ <1 eradicate disease from community given that 


the disease free equilibrium is globally asymptotically stable. 


3.5 Global Stability of Endemic Equilibrium Point of a model with intervention. 

In this section we prove the global stability of endemic equilibrium point E, of system (5) 
by using Lyapunov's direct method. Our Lyapunov function is constructed from suitable 
choice of logarithmic function. The global properties of endemic equilibrium point are 
studied by stating and proving the following theorem. 

Theorem 3.6 /f R, >1 then the unique endemic equilibrium E, of system (5) is globally 


asymptotically stable in the interior of Q. 
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Proof: We use approach of Korobeinikov [13] as it is used to most complicated 
compartmental epidemiological models, to construct the Lyapunov function from suitable 
choice of the following logarithmic function: 


W=} a(y,-y, In(y,)), 
where a, are properly chosen positive constants, y, is population of compartment i and 
y is the equilibrium level. We define the function 


W: {(s,v,L,i I h)eQ: S, v, l,i „i h>0}>R by: 


W (s,v,l,i,i,h)= A (s-s In(s))+A,(v—v'ln(v))+<A,(/-/'In(Z)) 
+A, (i, -i,In(i,)) +A, (i, -i,"In(i, ))+ A, (h—A'In(A)). 
The constants A,A,,---A. are non-negative in ( and W is Lyapunov function. The 
function W together with its constants A,A,,---A,>0 are chosen in such way that W is 
continuous and differentiable in a space C' and on the interior of ©, E, is global 
minimum of W on Q, and Wis, ie en ae ae h )=0. The time derivative of Lyapunov 


function W computed along the solutions of system (5) is: 


ds dv l \dl di, 
W'=A/1-— |—+A] 1-— 1-—|—+A,/1- 
ar; al ae al -Ja [iH ar 





(21) 
di, h’) dh 
+A, are A, (1-2 | 
At Endemic equilibrium point (EEP) we have: 
(I-p)A=(4+ B(i, +i, -8i -ôi )s" -v’, 
pi =(A+ 0-6, -ði Jv 
A+e= = Bs’ (i, +i, )+yph (i +i) + (di, +6,i, E), 
+6, += —=(I- nel’ +(1-$) oi, ) +64, +6, (22) 
4+0 +Ò, = Siar +0,i, , 
L, 
à= - (vi, + bei, —yB (i, +i ))+ di, +8. 


We re-write W” by using (22) as follows: 
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“fe = Mere ji, -ôi v -(2+0-4i, -ôi )v | 

saf EE] a +i,)+yBh(i, +i,)-—{ Bs’ (i +i )+ yBh (i +i, )+ (di) +44, E- (4), +6, yr} 
+A, (i a Ja- n)el +(1-) ai, -={(I- nel +(1-) i, +(6,i, +4, i, Ji — (6, + a 
afi 


er +(âi +43, ‘Vi, - (Gi, +â )i, bi, 
| + boi, = {oi + oi, —yB(i, +i, )+(5i, +6,i, )h +B (i, +i, )h -hôi + ain ph 


Simplification of (23) results to: 


Ai Ie ses a eae ad 


f 
W' = m ai = A 4 Ps, shi i h) . The function P(s,v,1,i i h) 
V 


balances the — hand side of (23). The function P(s,v,/,i,,i,,h) is non-positive following 


the approaches of McCluskey [15] and Mukandavire [17]. 
That is P<O for every s,v,l,i,i,h>0. Thus W'<0O for all s,v,l,i,i,,h>0O and zero when 


as Re 9% 9% 
s=s,v=v,l=I =0,i, =i, =0,i, =i, =0,h=h =0. Therefore the largest compact invariant 
set in 2 such that W'=0 is the singleton {£,' which is Endemic Equilibrium point of 


model (5). LaSalles's invariant principle [14] then implies that E, is globally asymptotically 











stable in the interior of the region 2 if R, >1and that completes our proof. 





4. NUMERICAL SIMULATIONS AND DISCUSSIONS 
In this section numerical simulation of normalized model (5) is carried out in order to 
illustrate the qualitative results by using available parameter values from existing 
literature as well as estimated ones. Unless otherwise stated parameter values appeared in 
Table 3 will be used during the simulation process. 


(23) 


20 / 23 
MLAY, LUBOOBI, KUZNETSOV, SHAHADA 


Table 3: Parameter values for normalized model (5). 


Symbol Value/range (yr` ) Source 

A 0.05 Estimated. 
p 2.58 Estimated. 
p 0.4 Estimated. 
0 0.1 Estimated. 
E 0.03 [8] 

1] 0.7 (0.7-0.95) [19] 

4 0.01923 (0.01-0.04) [2] 

Ô, 0.3 (0.07-0.365) [21] 

Ô, 0.2 (0.07-0.365) [21] 

D 0.6 Estimated. 
w 0.2 Estimated. 
Dv 0.3 Estimated. 
y 0.2 Estimated. 


4.1 Numerical Simulation of a model (5) in presence of intervention and TB. 
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Figure 3: Shows the dynamics of susceptible, vaccinated, latently infected, severely infected, 
mildly infected and treated population proportions in presence of interventions and TB 
with increasing time. 
Figure 3 shows dynamic behavior of susceptible, vaccinated, latently infected, severely 
infected, mildly infected and treated classes when R =1.8519. The plot is produced by 


MATLAB by using J =2.58;y =0.2;0 =0.1; € = 0.03;4 = 0.7; 2 = 0.05; 6, = 0.3; œ = 0.2; 
p =0.4;v =0.3;6, =0.2 6=0.6 as estimated parametric values and whose definitions are 


given in Table 2. Starting with initial values s(0)=0.60,v(0) =0.05,/(0) =0.1,i, (0) =0.1, 
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i, (0)=0.1 and h(0)=0.05, the system (5) attains the local asymptotic stability of endemic 
equilibrium point, E, =(s", v l iyi, h) = (0.4170, 0.1383, 0.3497, 0.0069, 0.0165,0.0716). In 


sh, sh, 
presence of interventions and TB, susceptible population proportion initially decreases to 
lower levels and later increases to its carrying capacity with time as shown in Figure 3. 
Vaccinated population proportion initially increases to higher levels and stabilizes as time 
increases. On the other hand both latently infected and treated population proportions 
increase to higher levels and gradually decreases to their carrying capacities. However both 
mildly and severely infected population proportions decreases to their lowest endemic 
levels. Again even with intervention, disease does not clear from community since effective 
reproduction number is R, =1.8519 >1. Classically this result supports the theorem of local 


stability of endemic equilibrium. 

4.2 Phase portraits illustrating dynamical behavior of population proportions at EEP. 

In this section phase portraits to illustrate the dynamics of the model (5) at endemic 
equilibrium point for susceptible class versus vaccinated, latently infected, severely 
infected, mildly infected, and treated classes are plotted by using parameter values 
indicated in Table 3. With different varying initial conditions, each solution curve in Figure 
4 tends to endemic equilibrium point E, presented in Section 4.1. Therefore we conclude 


that the system (5) is globally stable about endemic equilibrium point E, for the 
parameters displayed in Table 3. 


Vaccinated Proportion 











5 i i > i ' i i i 0 r r i i 
0.2 03 0.4 0.5 06 07 0.2 0.3 0.4 0.5 0.6 0.7 0.2 0.3 0.4 0.5 0.6 0.7 
Susceptible Susceptible Susceptible 


D 





Mildly Infected Proportion 











k i eee e k Susceptible 

Figure 4: Shows Phase plane portraits for dynamics of susceptible population proportion 
and (A) vaccinated (B) latently infected (C) severely infected (D) mildly infected (E) 
treated population proportions showing endemic equilibrium point with varying initial 
values as time increases. 
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5. CONCLUSION 
In this article, a continuous time deterministic Tuberculosis model with vaccination and 
treatment as intervention strategies has been formulated and the role of reinfection on 
transmission dynamics of TB is critically assessed. In presence of reinfection and multiple 
equilibria the backward bifurcation occurs at effective reproduction number R =1. In this 


scenario stable disease free equilibrium coexists with two endemic equilibria: smaller 
unstable endemic equilibrium (with small number of infected individuals) and larger stable 
endemic equilibrium (with large number of infected individuals) in the neghbourhood of 1 
when R, <1. This shows that even with classically reducing the threshold R, below one the 


disease still persist in the community. We suggest that reinfection is a real TB feature and 
an important aspect to consider when modeling the complex dynamics of TB. 
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